#define SEQAN_TRANSLATION_NO_PARALLEL
#include <iostream>
#include <fstream>
#include <seqan/basic.h>
#include <seqan/file.h>
#include <seqan/sequence.h>
#include <seqan/translation.h>
#include <seqan/align.h>
#include <seqan/graph_align.h>
#include <seqan/graph_msa.h>
#include <seqan/graph_types.h>
#include <seqan/store.h>
#include <seqan/seq_io.h>
using namespace seqan;

typedef String<String<Dna> > DnaVec;
typedef String<AminoAcid> TSequence;

bool IsG(Dna const & d)
{
    return d == (Dna)'G';
}

bool IsC(Dna const & d)
{
    return d == (Dna)'C';
}

void GCContent(String<Dna> dna)
{
    int GC = 0;
    for(unsigned int i = 0; i < length(dna); ++i)
    {
        if(IsG(dna[i]) || IsC(dna[i]))
            GC++;
    }
    std::cout << "GC-content is " << (double)GC/(double)length(dna) << std::endl;
}

void process(DnaVec dna)
{
    Align<TSequence> align;
    resize(rows(align), length(dna));
    typedef String<AminoAcid> TSequence;
    StringSet<TSequence> seq;

    for(unsigned i = 0; i < length(dna); ++i)
    {
        GCContent(dna[i]);
        String<Rna> rna(dna[i]);
        std::cout << rna << std::endl;
        String<AminoAcid> aSeqs;
        translate(aSeqs, dna[i]);
        std::cout << aSeqs << std::endl;
        assignSource(row(align, i), aSeqs);
        appendValue(seq, aSeqs);
    }
    globalMsaAlignment(align, Blosum62(-1, -1));
    std::cout << align << std::endl;
    Graph<Alignment<StringSet<TSequence, Dependent<> > > > aliG(seq);
    globalMsaAlignment(aliG, Blosum62(-1, -1));
    std::ofstream out("/out.phy");
    out << length(rows(align)) << " " << sumOfPairsScore(aliG, Blosum62(-1, -1)) << std::endl;
    for(unsigned i = 0; i < length(rows(align)); ++i)
    {
        out << row(align, i) << std::endl;
    }
}

int main(int argc, char** argv)
{
    DnaVec vec;
    for(int i = 0; i < argc - 1; ++i)
    {
        String<char> id;
        String<Dna> sequence;

        std::string filename = argv[i+1];
        std::ifstream in(filename.c_str());

        if(filename.substr(filename.length() - 3, 3) == ".gb")
        {
            read(in, sequence, Genbank());
        }
        else
        if(filename.substr(filename.length() - 6, 6) == ".fasta")
        {
            read(in, sequence, Fasta());
        }
        else
        {
            std::cout << "Wrong file" << std::endl;
            continue;
        }
        appendValue(vec, sequence);
        std::cout << sequence << std::endl;
    }
    process(vec);
    return 0;
}

